* Quarterly_county_sunsets.do
* 2014.02.07
* Last update 2014.11.05: updated to compute 2013 sunsets
* Creates a quarterly average sunset for each county in the QCEW data 

capture log close
set more off
timer clear 1
timer on 1
clear
set matsize 10000

local work "/DIRECTORY"

log using "`work'/logs/QCEW_county_sunsets.log", replace

* Reading data
use "`work'/data/qcew/QCEW_raw_quarterly.dta", clear

* Subsetting to 1 record per county
keep fips latitude longitude
duplicates drop

* Merging TZ offsets
merge 1:1 fips using "`work'/data/county_tz/county_time_zone.dta", keep(1 3) keepusing(tz_offset time_zone) nogenerate

* Expanding to 4 quarters - revised 2014.11.05 for inclusion of 2013 data
expand (365*24)+6 /* +6 adjusts for leap years */
bysort fips: gen long date = td(01jan1990) if _n==1
bysort fips: replace date = date[_n-1] + 1 if missing(date)
format date %td
gen tq = qofd(date)
format tq %tq

* Generating state FIPS
gen fips_st = substr(fips, 1, 2)

// Fixing historical time zones in indiana 
replace time_zone = "C" if (fips == "18027" | fips == "18037" | fips == "18083" | fips == "18101" | fips == "18125") & date < td(4nov2007)
replace time_zone = "C" if (fips == "18131") & date < td(11mar2007)
replace time_zone = "e" if fips_st == "18" & date < td(2apr2006) & time_zone == "E" /* this line modified from ATUS_proc to use state fips */

// Getting time zone offset right
gen observe_dst = 1
replace observe_dst = 0 if time_zone == "e" | time_zone == "h" | time_zone == "m"
gen dst = 0
rename tz_offset tz_offset_nodst
replace dst = 1 if observe_dst == 1 & date >= td(1apr1990) & date < td(28oct1990)
replace dst = 1 if observe_dst == 1 & date >= td(7apr1991) & date < td(27oct1991)
replace dst = 1 if observe_dst == 1 & date >= td(5apr1992) & date < td(25oct1992)
replace dst = 1 if observe_dst == 1 & date >= td(4apr1993) & date < td(31oct1993)
replace dst = 1 if observe_dst == 1 & date >= td(3apr1994) & date < td(30oct1994)
replace dst = 1 if observe_dst == 1 & date >= td(2apr1995) & date < td(29oct1995)
replace dst = 1 if observe_dst == 1 & date >= td(7apr1996) & date < td(27oct1996)
replace dst = 1 if observe_dst == 1 & date >= td(6apr1997) & date < td(26oct1997)
replace dst = 1 if observe_dst == 1 & date >= td(5apr1998) & date < td(25oct1998)
replace dst = 1 if observe_dst == 1 & date >= td(4apr1999) & date < td(31oct1999)
replace dst = 1 if observe_dst == 1 & date >= td(2apr2000) & date < td(29oct2000)
replace dst = 1 if observe_dst == 1 & date >= td(1apr2001) & date < td(28oct2001)
replace dst = 1 if observe_dst == 1 & date >= td(7apr2002) & date < td(27oct2002)
replace dst = 1 if observe_dst == 1 & date >= td(6apr2003) & date < td(26oct2003)
replace dst = 1 if observe_dst == 1 & date >= td(4apr2004) & date < td(31oct2004)
replace dst = 1 if observe_dst == 1 & date >= td(3apr2005) & date < td(30oct2005)
replace dst = 1 if observe_dst == 1 & date >= td(2apr2006) & date < td(29oct2006)
replace dst = 1 if observe_dst == 1 & date >= td(11mar2007) & date < td(4nov2007)
replace dst = 1 if observe_dst == 1 & date >= td(9mar2008) & date < td(2nov2008)
replace dst = 1 if observe_dst == 1 & date >= td(8mar2009) & date < td(1nov2009)
replace dst = 1 if observe_dst == 1 & date >= td(14mar2010) & date < td(7nov2010)
replace dst = 1 if observe_dst == 1 & date >= td(13mar2011) & date < td(6nov2011)
replace dst = 1 if observe_dst == 1 & date >= td(11mar2012) & date < td(4nov2012)
replace dst = 1 if observe_dst == 1 & date >= td(10mar2013) & date < td(3nov2013)
replace dst = 1 if observe_dst == 1 & date >= td(9mar2014) & date < td(2nov2014)
replace dst = 1 if observe_dst == 1 & date >= td(8mar2015) & date < td(1nov2015)
gen tz_offset = tz_offset_nodst + dst

drop if latitude == .
* Sunset times
solar_calculator latitude longitude date tz_offset
drop sun_ra sun_dec solar_noon
rename sunset_time sunset_time_avg
rename sunrise_time sunrise_time_avg
rename sunlight_duration sunlight_duration_avg

* Pretending no DST
replace tz_offset = tz_offset_nodst
solar_calculator latitude longitude date tz_offset
drop sun_ra sun_dec solar_noon sunrise_time sunlight_duration
rename sunset_time sunset_time_avg_nodst

* Mode of time zone for collapse
bysort fips tq: egen modal_time_zone=mode(time_zone)
encode modal_time_zone, gen(tz)

* Collapse to quarterly average
collapse (mean) sunset_time_avg sunset_time_avg_nodst sunrise_time_avg sunlight_duration_avg tz (min) tz_offset observe_dst, by(tq fips)

* Declare panel
encode fips, generate(panelID)
xtset panelID tq

* Saving
saveold "`work'/data/qcew/Quarterly_county_sunsets.dta", replace





timer off 1
timer list 1
capture log close



